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We propose a scheme for investigating the quantum dy- 
namics of interacting electron models by means of a time- 
dependent variational principle and spin coherent states of 
space lattice operators. We apply such a scheme to the one- 
dimensional Hubbard model, and solve the resulting equa- 
tions in different regimes. In particular we find that at low 
densities the dynamics is mapped into two coupled nonlinear 
Schrodinger equations, whereas near half-filling the model is 
described by two coupled Josephson-junction arrays. Focus- 
ing then to the case in which only the phases of the spin 
variables are dynamically active, we examine a number of dif- 
ferent solutions corresponding to the excitations of few macro- 
scopic modes. Based on fixed-point equations of the simpler 
among them, we show that the standard one-band ground- 
state phase space is found. 



I. INTRODUCTION 

Investigating quantum dynamics of strongly correlated 
many-body systems is a hard task since, even for ex- 
tremely simplified models, the interactions of the large 
number of degrees of freedom are usually affected by a 
nonlinear character. At the operational level this entails 
the impossibility of evaluating explicitly the action of the 
propagator known from the Schrodinger equation, that 
is the evolution |$) = exp[-itH/h]\^o) of a state |$o) 
governed by the Hamiltonian H. A standard way to re- 
duce such a difficulty to a more tractable form consists in 
recasting the purely quantum problem within an appro- 
priate coherent-states picture once the algebraic struc- 
ture characterizing H has been identified. This leads to 
represent the system evolution through the equations of 
motion issued from an effective classical Hamiltonian H 
expressed in terms of the coherent state parameters jl) . 

A systematic development of such an approach is 
provided by the time-dependent variational principle 
(TDVP) procedure [^. This amounts to constructing a 
trial macroscopic wave function \^) that contains time- 
dependent parameters whose evolution is derived so as 
to optimize the approximation of the quantum propa- 
gator action On this basis, using the generalized 
coherent-states to construct the trial state I^E") is quite 
advantageous in that the coherent state parameters nat- 
urally label l^*) and make explicit its dependence on the 
algebraic structure of H, namely on the operators de- 



scribing the microscopic physical processes therein. By 
making the phase that appears in j^*) coincide with the 
effective action, the Schrodinger equation turns out to be 
automatically satisfied when projected onto j^E'). 

In a recent paper Q such a scheme was specialized to 
the case of interacting electrons described by the Hub- 
bard Hamiltonian. There the coherent states entering 
I '5) were specific to the physical regimes (e.g., supercon- 
ducting, antiferromagnetic, etc.), the latter selecting case 
by case the appropriate approximate algebraic framework 
within the Hamiltonian dynamical algebra. 

The standpoint here adopted is instead to implement 
a unified TDVP treatment independent of the particular 
physical regime and provide a coherent state picture of 
electrons on the ambient lattice, whatever the model in- 
teraction actually is. Even though this approach is quite 
general, in the sequel we shall develop it for the Hubbard 
Hamiltonian. 

It is well known that the Hubbard Hamiltonian can 
be rewritten in terms of two coupled XX models of 1/2 
spin operators by means of the Jordan- Wigner transfor- 
mation. Such a transformation can be performed in any 
dimension as well as, in principle, for any electron Hamil- 
tonian, and leads quite naturally to a picture relying on 
spin coherent states (SCS) When this is used explic- 
itly within the TDVP framework, the resulting equations 
of motion are recognized to describe two coupled fluids, 
which dynamics we shall discuss. 

A basic trait of the spin description is that its semi- 
classical version is more reliable the more the spins are 
large [||. Since this feature is in general not realized 
when starting from quantum 1/2 spin operators, we shall 
look here, in particular, for solutions of the equations of 
motions corresponding to the macroscopic excitations of 
few system modes, in which case we expect to describe 
actual regimes for the Hubbard model itself. The prob- 
lem of mode requantization, naturally in order due to 
the expected quantum character of the low-temperature 
regime, is left to a successive analysis Q. 

The choice of j^*) as a direct product of single-site 
Bloch states, representing the only assumption for our 
construction, deserves some comments as to the expected 
reduction of the number of states in the Hilbert space 
that are actually available for the system dynamics. Such 
an effect usually occurs in a number of mean-field approx- 
imations like the standard Hartree-Fock (HP) in which 
the dominating features of the system are accounted for 
in an explicit way thanks to an extreme reduction of the 
states accessible to the system. 
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In this respect, using coherent states relative to the op- 
erators of H defined in the ambient lattice is by construc- 
tion less restrictive than using a subset of states tailored 
for a specific regime. The advantage coming from this 
choice is manifold. First, the structure of \^) is however 
able to produce an effective hamiltonian Ti that inher- 
its both the nonlocal and the nonlinear character of H , 
contrary to the Hartree-Fock (HF) scheme, in which H 
reduces to a sum of single-site linear Hamiltonians. In 
passing, we notice that in many cases 7i exhibits a form 
that is endowed with the same complexity oi H . In fact, 
the nontrivial form of Ti reflects the basic character of 
the TDVP method that singles out \^) variationally as 
the best solution to the original Schrodinger equation [|j , 
whereas within the HF approximation what is solved is 
a different Schrodinger equation, involving just the lin- 
earized Hamiltonian. 

Second, as a consequence of the above feature, also the 
propagation of any initial state is sustained by the full 
Hamiltonian, rather than by its linearized HF version. 
Indeed, it is easily shown that the latter entails quan- 
tum states whose time evolution is periodic, while the 
TDVP dynamics is endowed with a much richer struc- 
ture. In particular, the dynamics of the expectation val- 
ues of spin operators (our dynamical variables) is consis- 
tently reproduced, whereas, when turning to expectation 
values of products of spin operators, the description ob- 
tained does not differ substantially from the one that can 
be achieved within the random-phase approximation. 

The Jordan- Wigner transformation mentioned above 
amounts to rewriting the electron annihilation opera- 
tors Cj^ri, with rj ^tiii terms of Pauli spin matrices 
^a,j,Ta,j, with a = 1,2,3, which locally form two (com- 
muting) su{2) algebras. For the Hubbard model, it turns 
out that in dimension D > 1 the possible transformed 
Hamiltonians differ from each other due to a certain ex- 
ponential factor in front of the hopping term, which form 
in fact depends on the ordering chosen for labeling the 
lattice sites. This problem has been already investigated 
in the literature and in the present paper we shall 
limit our discussion to the one-dimensional (ID) case. 
Explicitly, 



Cj,T = Pjif^aWj 1 = PL(CT3)Pj(r3)r" 



(1) 



where Pj(i/3) = II^<j <^3.ii — '7,t, from which the ex- 
pressions for cl ^ are straightforwardly derived. Here L 



is the number of lattice sites, 



,+ 



' W2 7 J with h' = 



a,T. Remarkably, this transformation maps fermions, 
which anticommute on different sites, into spins, which 
commute on different sites, i.e. [ca.j, o'b^i] — for £ ^ j. 

Once Eqs. (|l|) are inserted into the Hubbard Hamilto- 
nian, the latter becomes 



H = ^[[7(73,, T,^, - Tia+a-^, + t+tJ^, + H.c.)] , (2) 
when periodic boundary conditions are considered, and 



an odd number of holes iV^ (ry =t,i) on both a sub- 
lattices is assumed 0, otherwise boundary terms (cor- 
responding to i = L) in the hopping contribution de- 
pending on T have to be rewritten as (e'^(^~^T ^cr^cr^ + 
g^(i-Jvf )^+^- jj^ y ^Yi^ gx^j.j^ ieYiaa that 

take advantage of conserved quantities such as the to- 
tal electron number and the magnetization have been ig- 
nored. 

In the next section, based on the spin-coherent- 
state picture, we shall implement the TDVP procedure 
whereby one can derive from (^) the effective Hamil- 
tonian and the related motion equations. In Sec. HI, 
upon recognizing the two-fluid structure of the resulting 
model, we shall solve explicitly the motion equations of 
each fluid within a phase-locking approximation, and ev- 
idenciate how the Coulomb interaction drives the system 
to a transition (apparently related to the metal-insulator 
one) in which also the phases of the two fluids become 
strongly locked. Tunneling phenomena between the two 
fluids are also discussed. In Sec. IV we specialize to 
the study of solutions exhibiting a pure phase dynam- 
ics, and stress the aspect concerning the macroscopicity 
of the excited degrees of freedom. In Sec. V we show 
that the ground-state phase space known from standard 
mean-field treatments can be obtained within our scheme 
by analyzing the fixed points of a very simple collective 
phase solution, corresponding, in fact, to describe the 
whole lattice as a sum of two-site clusters. Finally Sec. 
VI is devoted to give some conclusions. 



II. COHERENT-STATES PICTURE 

Approaching interacting spin systems within a semi- 
classical limit has been deeply investigated. In partic- 
ular, it is well understood that a consistent description 
can be obtained ||l| by projecting the Hamiltonian onto 
a basis of SCS. In this case, an exact result obtained by 
Lieb Q shows that the projected Hamiltonian reproduces 
the behavior of the original one the more the spins are 
large, and in any case it gives upper and lower bounds 
to the ground-state energy of the quantum Hamiltonian 
(the exact value being recovered for infinitely large spins) . 
One-half SCS are given by 

\ri) ^ c{i^y>'- 



1/2), 



(3) 



where the maximum weight vector | — 1/2) belongs to the 
J3 spectrum (J3I ± 1/2) = ±1/2| ± 1/2)) and fulfills the 
conditon J-|0) — 0, J_ (J+ = (J_)+) representing the 
lowering (raising) operator. Also, defining the normaliza- 
tion factor as C{ri) = 1/^1 -f jryp ensures the condition 
{ri\Tf) = 1. The expectation values of generators J3, J± 



S3 = (J3) = 



- 1 



2(1 + 



(1 + hp) 



(4) 



(5) 
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obtained by means of definition (||) ((•) = {r]\ • j?])), 
clearly exhibit their semiclassical character when consid- 
ering the fact that ^3 , S± satisfy the equation 8^ + 82 + 
Sf = 1/4 ((5*+ = Si + 182)), namely the same sphere 
equation fulfilled by the classic counterpart of the spin 
( Ji, J2, J3) ( J+ = Ji + i J2). In passing we notice that 
the spin variables, assuming limited values, keep track of 
the fermionic nature of the underlying system. 

The set-up just developed can be readily extended to 
the interacting spins of H. Assigning at each site a pair 
of SCS \aj), \(3j) relative to the above a-spin and r-spin, 
respectively, allows one to implement the TDVP proce- 
dure that is essentially based on constructing a macro- 
scopic trial wave function accounting for the microscopic 
processes of the system. The simplest choice for a spin 
model is realized through the state 



iA 



(6) 



where la) 



pectation values A* 



iSij{\aj) ® that provides the ex- 

= (B* = (*|t+|*)) and 

Asj = («'|cr3,il*) {Bsj = («'|t3j|*)) of a spins (r spins). 
The description of the microscopic dynamical activity in 
terms of such semiclassical variables (actually they cor- 
respond to an ensemble of classical spins) is achieved by 
showing that they obey a set of Hamiltonian equations 
standardly derived from imposing l^t^) to obey the weaker 
version of the Schrodinger equation {'i'Kihdt— H)\'i>) — 0, 
the latter requirement leading as well to interpret 8 in 
Eq. as the effective action. The explicit form of 
TDVP Hamiltonian generating such hamiltonian equa- 
tions turns out to be 

(H) = {P\(^{a\H\a)(g>\P) , 

while the Poisson Brackets obeyed by the spin ensemble 
variables implicitly follow from the equations of motion 
themselves. 

Hubbard Hamiltonian (0) in one-dimension, when pro- 
jected onto the trial state^a) (E" |/3), becomes 

(H) ^N,^ + ^{As + B,) + uY^ A3, Bs, + Ht , (7) 



where A^ = A^j , B^ = B^j and the hopping term 
TYt, which reads 

Ht - -r^(A*A,+i + B*B,+i + B.C.) , 



is nothing but the sum of two classical XX models. The 
Hamiltonian equations generated by the TDVP proce- 
dure are given by 

iA, = [-5a + UB3j)A, + 2TA3j{Aj^i + A,-i) , (8) 



iB, = {-6b + UA3j)Bj + 2TB3,{Bj^i + B,_i) , (9) 



3i 



T[A* (A,^, + ^,_i) - A, {A* + AU)] , (10) 



iBs, = -T[B* (5,^1 + i3,_i) - B, {Bl, + B*_,}] , (11 

where 6^ = fiA ~^ U /2, Sb = Hb — U /2, once the Hamil- 
tonian (H) is rewritten in the form 



H= {H) + ^.A Xa + A^B Xb 



(12) 



containing the constraints xa ^ o-a ~ A3, xb = tb — B3 
with Lagrange multipliers fiA , I^b ■ The Poisson brackets 
implicitly entailed by Eqs. (p|)-(fl]) are given by 



{c*,c,} = 2C3j/ih, {C3j,c*} = q/ih 

with C = A,B, and exhibit the structure of a (classical) 
angular momentum algebra. Also, they state that A3, 
B3 , related to the total number of spin- up and spin-down 
electrons by the formulas 



j 



) - A3 + N,/2 , ( 



j 



B3 + NJ2, 



respectively, where rijo- — c^^Cjc, (cr [), are constants 
of motion since {^3,?^} = = {^3,'^}. It is thus natu- 
ral investigating spin dynamics when A3, B3 are assumed 
to have fixed values va , ^b by inserting such informations 
via the constraints xa = = xb- 

The conservation, for each j, of the Casimir functions 
Caj = Alj + \Aj\'^ and Cbj = + \Bj\'^ is preserved 
as well. On the contrary, the total magnetization vector 



M = (M, 
with M+ 



^,My,M^) 



Mj (where M^, + iMy = M+ 



E,(4'l<T-|vE' 



^3 ^' '~ 3 '3 I ^ ' "3 1' 

conserved but only its z component — \ Ej {^\ 1^3 j ^ 

In addition, we also notice 



^jA*Bj) is no longer 



Ti:3)\^) = 2 Z^jA^3j - 033) 

that the usual particle-hole symmetry of the quantum 
Hamiltonian survives at the semiclassical level, and it is 
implemented by the particle-hole transformation A3j — > 
-A3j and B3J -Bsj. 

Two remarks are now in order. First, due to the choice 
of macroscopic wave function (|^), Hamiltonian (^, and 
Eqs. (H)-(ll ) mantain the same structure of Hamiltonian 
(^ and of the ensuing Heisenberg equations for the quan- 
tum spin variables, respectively, which feature is nontriv- 
ial §. 

Moreover, we notice that, when moving from the lat- 
tice description to the continuum limit |^ {Cj — > C{x) — 
|C(2:)|e^^(^), X & R, C = A,B), the resulting equa- 
tions can be interpreted as two nonlinear Schrodinger 
equations (NLSE) for the order parameter fields A{x), 
B[x). A part from the nonlinearity issued from C3J = 
±-\/l/4 — |Cjp that is capable of producing the stan- 
dard quartic term \Cj\'^ for |Cjp << 1/4, a further con- 
tribution in this sense comes from the Coulomb terms 
UA^jB^j. The standard reduction of the nonlinear 
Schrodinger equation to the continuity and the Bernoulli 
equation B governing the dynamics of the densitylike 
field |C(x)f and the phase field 0{x), respectively, sug- 
gests that Eqs. (p|)-(pT|) can be seen as describing the 
dynamics of a coupled two-fluid lattice model. 



5E,(^3 



B3 
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III. TWO-FLUID DYNAMICS 

The two-fluid structure of Eqs. has been rec- 

ognized by reducing them to the standard form (cubic 
NLSE) thanks to the assumption |Cjp << 1/4, namely, 
considering low-density fluids. In this regime the usual 
hydrodynamic picture is made far more complicated by 
the presence of A^j, B^iin front of the off-site T terms 
in Eq. (||), and Eq. (||). In fact such factors, in ad- 
dition to the usual Laplacian-like terms of the (lattice) 
Schrodinger equation characterized by A^j, B^j ~ —1/2, 
allow for the occurrence of configurations where the T 
terms exhibit anomalous signs (Asj, B^j > 0) through 
extended regions of the lattice. The investigations of the 
corresponding dynamics is deferred to a future study. 

A regime exhibiting, in a sense, an opposite character 
(|Cjp ~ 1/4 ^ Csj ~ 0) will be examined in the present 
section. The two-fluid structure still characterizes the 
motion equations even if the dynamics mainly concerns 
the phase variables, the densitylike variables |CjP be- 
ing now essentially constant. It is worth noting as well 
how such a regime (characterized by a Bernoulli-like dy- 
namics) is nothing but that the quantum phase regime 
naturally emerging from the XX model form of Tix for 
\Cj \ = const. In fact, by setting first 

Aj=Rjex.piaj , Bj—Sje^qiiPj , (13) 

where i?| = 1/4 - Aj^, S'j = 1/4 - B^j, consistently 
equipped with the standard canonical commutation rela- 
tions {ttf, ^3,} = 5ej/iH = {Pi, B^j}, and recasting then 
Eqs. (§)-(|ll|) in the action-angle variable version con- 
tained in the Appendix, one is able to work out the two 
linear second-order equations, 

dj = 4T'^[w{l3j+i - 2/3j + /3j-_i) + (a^+i - 2aj + a^^i)] , 

(14) 

fij = AT^[w{aj+i - 2a, + a^-i) + (/3j+i - 2(3, + /3,_i)] , 

(15) 

with w = U/AT, under the assumptions lA^j], \ B^\ « 
1/2, (a,+i - a,) « « (/3,+i - (3,). Eqs. (0-^) 
describe dynamics of first order quantities and exhibit the 
lagrangian structure typical of two classical planar XX 
models nontrivially phase coupled for any nonvanishing 
U ^0. 

Remarkably Eqs. ( [T^ ) and (|lj) can be decoupled (and 
solved) upon defining 9j = a, + (3j , {pj ~ a, — (3j . In this 
case they become 

- 4T2(1 + u.)(0j+i - 26, + 0,_i) , 
(p, = 4T^(1 - w){lp,+i - 2lp, + Lpj^i) , (16) 

whose solution can be easily worked out in terms of 
Fourier modes. In particular, let us notice that the pa- 
rameter w plays a relevant role, in that it drives the 



ip dynamics of the system from an oscillatory regime 
[w < 1) to a damped one {w > 1), whereas the 9 
dynamics remains purely oscillatory. This is explicit 
when considering any single mode solution of the form 
ipj{t; q) — cos(Ag< -I- i^j) and the ensuing dispersion rela- 
tion 

Xl = 16T2(1 - w) si-a^inq/L) . (17) 

In terms of the original phases a, and (3j this implies a 
phase- locking phenomenon for w > 1 (U > 4r), which is 
physically quite natural the more the on-site Coulomb 
repulsion becomes large. Having in mind the metal- 
insulator transition typical of the Hubbard model, which 
takes place at analogous values of U, we can argue that 
the change in the dynamical behavior parametrized by w 
might bear memory of such transition. 

It is worth noting that, again to the first order, Eqs. 
( |lO| ) and ( |Tl| ) for A^j, B^j reduce to 

is, - -(T/2)(a,+i - 2a, + a,.,) , (18) 

i?3, = -(r/2)(/3,+i - 2f3, + /3,_i) (19) 

which, despite the approximation introduced, still shows 
a nontrivial time dependence of A^j , B^j . The compari- 
son of the above equations with those descri bing the tun- 
neling phenomena of Josephson junctions |11| is quite 
natural, coming from the fact the same equations can 
be obtained, in the same linearized form, when con- 
sidering the Josephson-junction array Hamiltonian that 
can be represented in the sim plifi ed form by Hjj — 
Ej C'lj - 5 cos (7j+i - 7j) This is confirmed 

as well by Eqs. and (^ of the Appendix which, 

within the present approximation {Rj,Sj ~ 1/2), repro- 
duce exactly the equation C^j = {C^j, Hjj} for the on- 
site charges Csj . The special trait characterizing H is the 
quadratic term A^-jB^-j that generates a coupled phase 
dynamics via Eqs. ( p^ and (^5|), namely a linearized 
system of two [/-coupled arrays. Also, this suggests to 
define here a quantity that describes the net local current 
between the two arrays. If we let Aj and Bj play the role 
of the Josephson wave functions, and A^,, B^, as on-site 
charges, such current turns out to satisfy the equation 

T 

--^i^Pi+i-'^fj + Vj-i), (20) 

where Ij = A^j — B^j . Hence the tunneling phenomenon 
keeps track itself of the dependence on w, vanishing in 
the strong Coulomb repulsion regime {U > AT). 

IV. PHASE DYNAMICS 

Apart from the case related to Eqs. ( p^ and (p^, in 
the present paper we shall investigate solutions of Eqs. 
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such that only the phases play a relevant dy- 
namical role, A^j and B^j being constant in time. If, on 
the one hand, the dynamical situations in which A^j , B^j 
are involved exhibit a complex behavior and their inves- 
tigation goes beyond the purposes of the present paper, 
on the other hand, considering only aj, (3j as dynami- 
cally active still entails situations that are far from being 
trivial and facilitates the recognition of the topological 
features that possibly characterize the solutions. 

Hamiltonian describes the dynamics of interacting 
classical angular momenta. The latter exhibits solutions 
that consistently match the semiclassical nature of the 
present approach the more, by appropriately changing 
the basis of canonical coordinates, one identifies some 
new variables that could assume macroscopically large 
values and exhaustively account for the system dynamics 
1^^ . In general, for a given dynamical system, the ex- 
citations corresponding to the proper dynamical modes 
(if any) provide both the simplest and natural way to 
construct macroscopic semiclassical solutions. Unfortu- 
nately the identification of proper modes is equivalent to 
making explicit solution of the Hamiltonian equations, 
which in our case are highly nonlinear. Nevertheless, 
based on the usual Fourier modes picture, where 



=L-i/2^exp(ifcj)C'fc , 

k=l 

with k = i'Kk/L^ C — A,B, one may wonder whether 
there exists any integrable case corresponding to asso- 
ciate the macroscopically large number of spin degrees of 
freedom with a finite number of excited Fourier modes. 
It turns out that this is the case at least for two classes 
of solutions. 



A. Vortex dynamics 

First, it is easily verified that the case corresponding 
to two single excited Fourier modes p and q, one for each 
fluid, i.e., Ap = L^/'^Ra, Ak = 0,k p, and Bq = 
L^I'^Rb, Bk^O,k^q, is solution of Eqs. with 



Ep,g = U{A +1/2){B + 1/2) - 2T[R\ cos p+Rl cos q] . 

(25) 



Aj{t) = RAexp{i[{jp - ujA{p)t + (Pa)]} , 
Bj{t) = RBe-x.p{i[jq - ujB{q)t + (f>B)]} , 



(21) 
(22) 



where Rc = -^{-C^ with C ^ A,B, and A = Aj,/L, 

B = B3/L, (pA, (t>B are arbitrary phases accounting for 
the U{1) symmetry of dynamical equations and 



l^Aip) 



UB) +ATAcosp , 



UBiq) = i-SB + UA)+ATBcosq . 



(23) 



(24) 



The corresponding energy per site is straightforwardly 
obtained as 



The main feature of solutions ( |21| ) and (^2|) is their topo- 
logical character encoded by the winding numbers p and 
q. Notice that we have assumed periodic boundary con- 
ditions providing our ID lattice with the topology of the 
circle, and Aj, Bj can be regarded as order parameters 
covering two configuration spaces. Within this pic- 
ture the indices p and q account for the number of times 
Aj and Bj cover their configurations spaces while j goes 
from to L. Indeed such configurations are nothing but 
ID vortex excitations once the phases of the order pa- 
rameters are identified with the potential functions of 
two coupled fluids. Here the coupling is fully contained 
in the frequencies u!a{p) and ujB{q)- 

Interestingly, it is possible to evaluate explicitly corre- 
lation functions for solutions (|l|) and (H). Their physi- 
cal meaning is better understood when writing them for 
the original fermionic system. In this case, two-site cor- 
relations within a single fluid (e. g., the one with up 
spins), read 

(c]^c,t + H.c.) = 2{2A)\'-'\-^Ra cos[p(j - ^)] , ^ ^ j , 

(26) 

whereas for sites belonging to the two different fluids are 
(cj^QX + H.C.) = 2{2A)''-^{2BY-^RaRb 
X cos{jp -eq+ [ujB{q) - i^A{p)]t + (0A - 0b)} • (27) 

with j 7^ £. As expected, in both cases long-range or- 
der does not emerge since 2\A\, 2\B\ are smaller than 
one in any nontrivial case. However, two remarkable fea- 
tures emerge. First, they manifestly keep track of the 
topological character of the solution through the winding 
numbers p and q. Second, but more important, the two- 
fluid correlation function also exhibits a time-dependent 
behavior, whenever the density of the two fluids or the 
topological charges are different. This last feature should 
be viable to experimental observation. 



B. Staggered dynamics 

The general class of solutions characterized by the 
phase dynamics is obtained when B^j, A^^ are assumed 
to be assigned. In this case Eqs. (^)-(|ll]) reduce to a lin- 
ear system of equations for the variables Aj^s, and Bj's 
where proper modes coincide with the eigenvalues of a 
certain secular equation. In fact, one should recall that 
assigning B^j, A^j and thereby reconstructing \Bj\, \Aj\, 
leaves the possibility to satisfy the eigenvalue problem by 
exploiting just the phases of Bj and Aj. 

For A^j and B^j constant in time, Eqs. (|o|) and ( [Tl| ) 
are conveniently rewritten (see the Appendix) in terms 
of action-angle-like variables defined in (O), as 
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Rj+1 sin (a 



•i+i 



Rj-i sin(aj_i — Qfj) = , (28) 



Sj+i sin(/3,+i - + 5,_i sin(/3,_i - = . (29) 

The general solution is not known. Of course a simple 
solvable case is obtained by assuming both Rj and Sj 
constant and independent of j. This leads to the vor- 
tex case discussed in the previous subsection. A further 
solution exhibiting an interesting dynamics is obtained 
by noticing that Rj+i, Rj^i, can be factored out from 
the above conditions upon assuming that R21 = Re and 
i?2«+i = -Ro, V-^, with Re, Ro fixed constants. The 
same assumptions can be implemented on S'j+i, Sj_-i, 
so that when they are inserted in Eqs. ( p8| ) and (29), 
these turn out to depend only on the difference 7j+i — 7j, 
with 7 = a, p. The latter has two possible values satis- 
fying the equations, 7* or tt — 7* for each j, with 7, 
time-dependent function. Then Eqs. (||) and (^) can be 
solved explicitly, when rewriting them in the action-angle 
form of the Appendix. In fact, it turns out that a consis- 
tent solution is achieved provided 72^+1 — 72 j = 7*, and 
l2j — l2j-i = TT — 7,, for each j, which entails 



l2j + l 



: jTT -I- 71 , 72j = (j - 1)tt + 72 



(30) 



71 and 72 are time-dependent functions responsible for 
the system's phase dynamics as solutions of the corre- 
sponding equations given in Eqs. (58), (|60|). For instance 
in the case 7 = a they read 



ai = {5A-UB^o)t + ai{Q) 
a2 = {6A-UB^E)t + a2{Q) 



(31) 



while the analogue for (32 is easily derived. Interest- 
ingly, the time-dependent part of the phases keeps track 
of the coupling between the two fluids for any nonvan- 
ishing value of the Coulomb repulsion U . Again, such a 
feature should be viable for experimental observation. 

Apart from the initial conditions 71 (0), 72(0), the solu- 
tion ( |30| ) , (^ij) clearly exhibits a staggering in the phases 
both on the even and on the odd sublattices. Making 
such a solution consistent with periodic boundary condi- 
tions constrains the length of the lattice L to be L = 4p, 
p G N. Once more this feature can be related to the 
macroscopic excitation of some Fourier modes (two for 
each fluid). Explicitly for C = A 



A 



L/4 



A 



3L/4 



REe"'^^^^ + iRoe'"'^°'> 



(32) 



and Ak — for k ^ p, 3p, the analogue holding as well 
iorC = B,(j) = (3. 

The minimum energy per site Eg of the above stag- 
gered solution -to be compared with successive results 
for different phases- is found to be Es = U{v + |i^|)/4. 



It is important to observe how the independence of Eg 
from T (to be interpreted as the absence of a net global 
current) follows from the fact that the contributions to 
the hopping term coming from subsequent lattice bonds, 
let us say (j, j -I- 1) and (j -I- 1, j -I- 2), cancel each other. 
At the microscopic level, however, the hopping terms ac- 
tually contribute in terms of local currents (these are es- 
sentially given by 72^+1 - 72^ = 7*. l2j - l2j-i = tt - 7*) 
with opposite sign. 



C. Many-sublattices solution 



Further solutions to Eqs. (^-(11) that corresponds 
to the excitation of a finite number of Fourier modes 
(endowed with a macroscopic character) can be recovered 
by partitioning first the lattice A into n ^ L/q sublattices 
Aa of q sites {q divisor of L), and introducing then the 
collective variables 



q-l 



9-1 



^a=Y^ Ain+a , = ^ A 



3{£n+a) 7 



(33) 



i=0 



with ^ e (0, g - 1), a e (1, n). Here Aa, A*, and A^a stiff 
fulffU the commutation relations of a (classical) algebra 
su{2). It turns out that Eqs. (p|)-(pT|) can be rewritten 
in terms of the above collective variables provided fur- 
ther assumptions are stated. These are A^j ~ A^a/q, 
Aj = Aa/q with j g Aa. When this is the case, dy- 
namical equations reduce to a set of 4,L/q equations now 
written in terms of ^3^, Aa, B^a, Ba exhibiting the same 
structure. In the Fourier transformed space this amounts 
to the excitations of n modes, i.e.. 



1 



n 
a=l 



'Aa , At 



0. 



(34) 



(k = 2T:k/N) for k — mq, I ^ mq (0 < m < n), respec- 
tively. Solutions within this class are now obtained by 
solving the remaining An equations, which preserve the 
same complex structure of the original ones. 

For the simplest case n = 2 (n = 1 being a subclass of 
vortexlike solutions) the dynamical equations are repre- 
sented by 



iAi 


- {-Sa - 


f UB3i)Ai - 


f 4TA31A2 


(35) 


iBi 


= {-Sb - 


hUA:,i)Bi 


f ATB31B2 


(36) 


iM 


= {-Sa ^ 


f UB32)A2 - 


f ATA32A1 


(37) 


iB2 = 


= {-5b 4 


- UA32)B2 4 


- 4TB32B1 , 


(38) 



together with those for A* and B* . Correspondingly 
Hamiltonian (|^) takes the form 



N. 



H2 = ^{U/2~ [5c{C3^ + C32)-vc] + 



C=A,B 



6 



f/(A3iB3i + ^32^32) - 2t{AiA; + Bib; + c.c.)} . (39) 

As the number of first integrals of motions is 3 {H, A3, 
and -B3), whereas the equations are now 8, this case is 
nonintegrable. However, being interested in phase dy- 
namics in which case A^j and B^j are constants for 
each j, the solution to Eqs. (^-(^ can be worked 
out explicitly. The latter is characterized by collective 
frequencies Xa, Ab for the Aj's and Bj's of the form 
Cj = Cj{0)exp{iXct) (C = A,B, j = 1,2) which are 
independent from each other. 

It is important to notice how the case presently studied 
differs from the staggered solutions described above since 
Cj+2 = Cj is not contained in Eqs. (|2^)-(p9|), namely 

0. When Cj{t) are inserted in 



one is able to recast them in the form 



U{B3l - B32) = 4r ( ^32^ - ^31^ 



U{A3i - A32) = 4r ( B32^ - Bsi^ 

±>2 -Dl 



25a - - Ui^B + 4r ( ^32 4^ + A3i^ 



(40) 



(41) 



(42) 



fixed points of the equations of motion. Many other 
(mean-field-like) approaches are focused on the study of 
equilibrium and especially ground states of Hamiltonian 
(fr). For instance, from the Hartee-Fock approximation 
p3[ it is known that the T = phase space contains 
an antiferromagnetic, a ferromagnetic, and a paramag- 
netic phase for U > 0. In this section we shall see that 
a similar description of the ground-state phase can be 
already obtained by studying fixed points of the simple 
two-sublattice solution, the latter being given by Eqs. 
(g^)-(||) for Ac = 0. 

In particular, as Eqs. ( ^2| ) and (^3| ) just foi Sa, Sb, we 
search for the solutions of Eqs. ([4Cl|) and (|4l|), in which 
the unknowns are two. It is convenient to introduce the 
pair of new variables a = A31A32, b = _B3ii332, in terms 
of which the above equations reduce to the pair of fourth 
order, coupled equations 



(-1 



Ab) 



(1 



4a) 2 - Ai^j^ 



{v\ - 4a) = g- 



{\^4bf(yl-4b) 
(1 + 46)2 



44 



(45) 



(46) 



with g = AT/U. The two equations (^) and (|^ can be 
recast into a single eight-degree equation for the variable 
a, 



25b -Xb = UVA + 4T B: 



>32 



Bi 
B2 



B31 



B2 
Bi 



(43) 



where Cj (C^j) stay for initial conditions Cj(0) ([C3j(0)], 
and the constant of motion 



VA = A31 + A32 , VB = B31 + B; 



32 



(44) 



are input data, whereas 5a, 5b, A^j, B^j (consistently 
with vatI^b = const) are the unknown parameters to be 
fixed. 

It is worth noticing that Eqs. (^) and (^TJ) turn out 
to be completely independent from Xa, Xb while in Eqs. 
( ^ ) and ( p3| ) A a and Ab can be incorporated inside 5 a 
and 5b by redefining them as Ac = 5c — Xc/2, C = A,B. 
At the operative level this fact allows one to reconstruct 
the solution of (40)-(^3|) for Ac ^ from the case Ac = 0, 
which by the way identifies the fixed points of Eqs. (^5|)- 
( p8| ) . The investigation of such points is deepened in the 
next section. 



V. FIXED POINTS OF TWO-SUBLATTICE 
SOLUTION 

The present approach is able to give a (simplified) de- 
scription of the system dynamics, with a number of in- 
teresting features, as we have seen in the previous sec- 
tion. Nevertheless, as a secondary effect, it also gives the 
system equilibrium states, which coincide in fact with 



[Z\l-g')-4^l][{l + ul){Z'-4,.l)- 



B) 

2 I r72 



g'Z\v\ + 1 - Z)Y = AviiZ' - \y\Y , 



(47) 



where Z = 1 + 4a, and the factor v\ — 4a, which provides 
an independent solution, has been factored out [see (i) 
below] . The variable b is then easily worked out from 
Eq.®. 

First, let us notice that the independent solution v\ — 
4a = imphcs lyg—Ab = leading to a = ^'1/4, b — v%/A- 
As 1/ = VA + t^B = 2{sa^/a + Sb^/b) with Sa — ±1, 
Sf, = ±1, this solution implies A31 — A32 = i^a/2, 
B31 — B32 — vb/2 (two-sublattice solutions with fer- 
romagneticlike order on each sublattice). It has energy 



^ H 



u 
2" 



1 + V + VAVB + |(i^l + 4 - 2) 



(48) 



which matches the one of vortex solution Ep q [see Eq. 
(p5|)] in the untwisted case p = g = 0. For fixed filling 
V = va + vb its minimum value depends on the actual 
value of p. When g < \ (i.e., > 4, this is reached either 
for j/^ = , = or for j^b = , j/^ = 0, in which case 
the energy becomes 



pf 



U 
2" 



1 



9 



9 2 

— V 

2 



(49) 



the solution describing ferromagnetism away from half- 
filling within a single cluster, in that the average magne- 
tization on the cluster M — [ua — vb)I2 coincides with 
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±1^/4. The subindex p in Epf is to remind us that the 
solution on the lattice, due to the arbitrary choice of the 
sign of M on each cluster, does not exhibit ferromagnetic 
order. 

On the contrary, for <? > 1 the minimum value of ex- 
pression (^ ) is reached when va = vb = v /2. Physically, 
it corresponds to a paramagnetic solution even within a 
single cluster, and has energy 



E-n — 



U 



(1 



5(1 



-) 



(50) 



Having in mind the phase diagram known from mean- 
field-like Hartree-Fock treatment of the Hubbard model, 
an antiferromagnetic solution is also expected, where the 
energy should be lower than both Ef and Ep near half- 
filling and for large U . This can be worked out as a solu- 
tion of Eq. ( ^7| ) when the magnetization is zero, naniely 
va = vb- In this case, it is easily realized that Eq. (|47| ) 
can be rewritten as the product of a second-order factor 
{g^ — 1)Z^ + Av^ (real for g < 1) and a sixth-order fac- 
tor which, in the range of parameters physically allowed, 
never provides real solutions. On the contrary, the van- 
ishing of the second-order factor in a, in fact, leads to 
an antiferromagnetic solution. This can be seen by first 
realizing that an analogous equation holds also for 6, so 
that finally Eq. (45) and (|6|) reduce to two second-order 
ones 



1 = 



1 = 



g^(l + 4a)^ 



4-1 



g^(i + 4b)^ 

(l + 4&)2 



(51) 



(52) 



which in order to consistently match va — >^b imply a = 
b, with 



1 



2Wa\ 



(53) 



Notice that, away from half- filling (which corresponds to 
I'A = = 0), the condition g < 1 follows from Eqs. 
( pT| ) and (^2|). Moreover, when calculating explicitly A^i 
(A32 = VA- ^32) and B31 (B32 ^ vb - B32) through 



formula (53) which read 



^31 = h^A ± 




1 - 2 



vb ± 




9' 



(54) 



(55) 



Wb\ < 



one singles out the further restriction {vaI 
\/ (1 — ,9)/(l -I- g). The apparent freedom in choosing the 
sign in Eqs. (ISJ) and (El) just corresponds to exchange 



the role of ^31 and ^32 (B31 and ^32). In fact, the 
physical solutions turn out to be just two, in that the 
condition a — b can be implemented in two different 
ways, namely ^31 — B^i, A32 — B32 (paramagnetic), 
and A31 = B32, A32 = S31 (antiferromagnetic). The 
energies corresponding to such solutions, 



U 



K = iTi^ + '^) + U{Al+Al2)-8TJ- - AU - - Ai 



1 



1 



^32 I 



2UA,,A,2-8TJ-^AU -~Ai 



^32 I 



differ only due to the U term, which is manifestly lower 
in the antiferromagnetic case. It turns out that the anti- 
ferromagnetic cluster energy 



E, 



U 



(56) 



is always lower than paramagnet ic case within t he do- 
main specified by {vaI — \i^b \ < \/(l — g)/(l + g). 

The successive comparison among Eqs. (EO), (49), and 



(p6| ) shows that the ground-state phase space for this 
two-sublattices solution (Fig. 1) exhibits a structure in 
qualitative agreement with many other theories, in par- 
ticular the one obtained in the low-density approxima- 
tion for the one-band model (Ref. jl^). Moving from 
half-filling, in which case a magnetic phase is found for 
U > 4T, the antiferromagnetic phase takes place at in- 
creasing values of U, and in any case for filling greater 
than one quarter. Indeed, by requiring that Eaf < Epf, 
the transition line to the antiferromagnetic phase is given 
by 



= - (g 
g ^ 



i + yr^ 



> 0. 



(57) 



For lower values of filling, the system is a nonmagnetic 
metal. Within such regime an extra transition emerges 
for g = 1 from a paramagnetic solution with ferromag- 
netic structure on each cluster (energy Epf), and a para- 
magnetic solution with no order even within the clusters 
(energy Ep). Apparently by increasing U the lattice be- 
gins to organize towards ferromagnetism. Let us notice 
that, consistently with the ID character of the model 
studied, both ferromagnetic and antiferromagnetic solu- 
tions exhibit only local order, in that the actual value of 
the magnetization on different two-site clusters is uncor- 
related. 

In the previous section we explicitly gave the energies 
corresponding to some simple solutions of the equation of 
motions exhibiting nontrivial dynamics. A natural ques- 
tion is then whether some of these solutions survive down 
to the ground state, or not. Interestingly, one can verify 
that, in fact, the staggered solution, with energy Eg, at 
half-filling turns out to be degenerate with the two-site 
antiferromagnetic solution described above, with energy 
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Eaf. Indeed, both of them in this case have a vanishing 
hopping term, in agreement with the expected insulating 
behavior of such regime, and in practice on the single 
cluster the two solutions coincide. However, the explicit 
solution of the equations of motion in the staggered case 
proves that at a dynamical level the only consistent way 
of moving from the fixed point is by means of the stag- 
gered choice of phases. 

VI. CONCLUSIONS 

The main object of this paper has been to develop an 
approach to the Hubbard model quantum dynamics that 
is not based on the particular physical regime under in- 
vestigation, on the one hand, and is capable of reformu- 
lating the model dynamics in a form more tractable than 
that relying on the direct diagonalization of the model 
Hamiltonian, on the other. Such requirements have been 
achieved by combining three ingredients, which are the 
representation of quantum dynamics within a coherent- 
state picture, the expression of the Hubbard Hamiltonian 
in terms of spin variable (^) issued from its fermionic 
standard form through the Jordan- Wigner transforma- 
tion (^, and the implementation of the TDVP method. 
The choice of the trial state (^) has generated Hamil- 
tonian (|^) (that is Ti. with the constraints xc = 0; 
C — A, B) whose dynamics is governed by Eqs. (^)- 
(|ll|), and accounts for the evolution of the spin operator 
expectation values. 

The resulting dynamical scenery has revealed both a 
rich structure -that corresponding to a pair of XX models 
coupled through the Coulomb term- and interesting links 
with other models. 

For \Aj\'^, \Bj\'^ ~ 0, one obtains a model of two cou- 
pled fluids at low density. In particular, in this limit, 
the dynamics has been recognized to have the form of 
two coupled lattice NLSE. A feature that is unusual for 
the standard NLSE comes from the dependence of the 
off-site terms in Eqs. (^) and on the signs of C^j, 
which allows for the fragmentation of the planar lattice 
in regions where either C^j ~ +1/2 (sites occupied by 
electrons of type C), or C^j ~ —1/2 (local depletion of 
electrons of type C). The latter case suggests the occur- 
rence of soliton-like behavior in correspondence to the 
negative sign of the off-site terms. 

The opposite regime \Bj\^ ~ 1/4 has been stud- 

ied in Sec. HI, where the equations of motion have re- 
vealed that the model actually describes two coupled 
Josephson-j unction arrays. In particular the condition 
|Cjp ~ 1/4 makes emerge a dynamics concerning essen- 
tially the phases aj, (3j that can be solved exactly af- 
ter reducing the equations of the Appendix to the linear 
system described by Eqs. ( p^ and (^5|). Its main fea- 
ture is certainly the macroscopic effect of phase locking 
[{aj~ Pj) 0] which is enacted when going from U < AT 
io U > 4T, and might be related to the metal-insulator 



transition exhibited by the Hubbard model. 

Pursuing the investigation of dynamical situations in 
which phases are active and = const has led 

us to recognize two other interesting results. First, a set 
of topological solutions has been obtained by considering 
uniform configurations Csj — C3/ C — A, B, which are 
nontrivial when one excludes the half-filling case. The 
phases aj and Pj are allowed to change as j is varied 
so as to give rise to a pair of vortexlike configurations 
labeled by two integers p and q for Aj and Bj (the fluid 
order parameters), respectively. Also, the time behavior 
exhibits a dependence on the electronic fillings as well 
as on the topological characters through the frequencies 
uJAip) and ujBiq)- 

A second class of solutions has been obtained instead 
when considering the solutions of Eqs. (||) and (j9|) ful- 
filling the constraints Csj = const at each site, and de- 
pending on a unique frequency. Despite the strong sim- 
plification thus introduced the complexity of the problem 
is still dramatic as shown by the dynamical constraints 
( p8| ) and (|2^). It is worth recalling that their implemen- 
tation corresponds to find first the eigenvalues of (|^) and 
(^ in which C^j are regarded as constant, assigned pa- 
rameters, and singling out then the subset of eigenvectors 
such that \Cj \ are compatible with the assigned Csj. The 
staggered solutions [see Eqs. (|3l| ) and (|3^)] represent the 
case where the avalanche of initial conditions is reduced 
to a set of four data namely the values of \Cj\ for the 
sublattices of both even and odd sites. 

Based on the polygonal symmetry of the spin equations 
of motion their number has been reduced by introducing 
the collective variables ( ^3|) in Sec. HI C. The first non- 
trivial case (but also the only one directly tractable in 
an analytic way) has been shown to correspond to the 
two-sublattice solutions {Cj — Cg, C^j = €31, I = i -I- 2, 
y?') . The analysis of the fixed points of Eqs. ( pq ) and 
( p8| ) allows one to reconstruct the set of configurations 
in which those corresponding to the minimum energy are 
implicitly contained as a consequence of the absence of 
dynamics. In Sect. V we specialized to the latter in order 
to obtain a zero-temperature phase space. Noticeably, 
we have seen that already such a simple two-sublattice 
solution contains all the qualitative features of similar 
diagrams obtained in many other theories. Hence, we ar- 
gue that the general solution of fixed-point equations, if 
available on finite lattices by means of numerical analysis, 
should exhibit a richer structure than the one obtained 
within standard mean-field schemes even for what con- 
cerns the zero-temperature phase-space. 

Further developments of the present work can be envis- 
aged along the following lines. As to the transformation 
(|l|) it is important to notice how its use has been possible 
because of the ID character of the system. In higher di- 
mensions, in fact, this transformation depends explicitly 
on the ID path employed to cover and thus enumerate ex- 
haustively the lattice sites. Such a dependence introduces 
in the hopping term of the Hamiltonian a site-dependent 
exponential phase factor, which does not prevent the im- 
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plementation of the approach developed here. Hence, in 
spite of the increased complexity thus introduced, a nat- 
ural extension of the present work is in the analysis of 
the 2D case dynamics. 

As a matter of fact, due to the large number of de- 
grees of freedom involved, the ID case itself is already 
not fully tractable via numerical investigations. In this 
respect, focusing on zero-dimensional systems is almost 
expected in order to have a dependable numerical de- 
scription. On the other hand, it is well known that the 
physics of such mesoscopic systems (e.g. quantum dots 
and Josephson junctions) is properly depicted in many 
cases by Hubbard- like Hamiltonians [Q, Further 
investigation of such systems within the scheme proposed 
here seems promising. 

A final point still deserves to be deepened, which is the 
requantization of the spin variables. Despite the obvious 
difficulty of such a task in general |^ , the dynamical sit- 
uations here investigated, involving the macroscopic ex- 
citation of few system modes, seems quite feasible to this 
end. 
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FIG. 1. The ground-state phase diagram of the one-band 
model for the two-site solution: d — 1/ is the electron doping 
{d — half-filling) on the two-site cluster, and 1/g = U/4T. 
Its structure is in qualitative agreement with the diagram of 
Ref. [10], p. 256. 
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where (j) indicates the set of the nearest- neighbor sites. 
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iTiath- 



APPENDIX 



After setting Aj — Rj expiaj, Bj = Sj expi(3j, where 
Rj = {1/4- Al/y^^, Sj EE {l/A-BijY^^ with the Poisson 
brackets {ae^A^j} = 6ij/ih, {PijB^j} = 6ij/ih, it is 
found 

R 

a J =6a- UB^j + 2TA3J cos(a, - a, ) , (58) 

A^j = 2TRj Y Rtsm{aj - a,), (59) 

=6b- UA^j + 2TB3, V 1^ cos(/3, -Pj), (60) 
ie(j) ^ 
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